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Computing  Block- Angular  Karmarkar  Projections  with 
Applications  to  Stochastic  Programming 
John  R.  Birge*and  Liqun  Qi  ** 


Abstract :  We  present  a  variant  of  Karmarkar’a,  algorithm  for  block  angular  structured  linear  programs, 


such  as  stochastic  linear  programs.  By  computing  the  projection  efficiently,  we  give  a  worst  case  bound  on 
the  order  of  the  running  time  that  can  be  an  order  of  magnitude  better  than  that  of  Karmarkar's  standard 
algorithm.  A  related  variant  is  applied  to  the  dual  program,  and  its  implications  for  very  large-scale  problems 


are  given. 
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1.  Introduction 


Block  angular  and  dual  block  angular  linear  programs  arise  in  a  variety  of  applications,  in  particular, 
stochastic  linear  programming.  Dual  block  angular  linear  programs  (see,  e.g.,  (10j,[21j,[29j)  have  the  form: 

T  ,  T 

mux  c£z0  +2Jk=icfc  *k 

subject  to  A0*o  =  b0,  .  . 

•dfc*o  =  4fc,A:  *  1, 

**  >  0,  4  =  0, . . . ,  N, 

where  the  sises  of  the  matrices  are  consistent  with  xk  6  R"* ,  k  =  0, . . . ,  N,  6*  €  #"*“ ,  k  =  0 .,N,  where 
mfc  <  Ufc,  k  —  0, . . . ,  N,  and  we  assume  Aq  and  Wk  have  full  row  ranks.  We  also  assume  that  the  blocks  are 
large  enough  that  no  <  52k=l  nk. 

Stochastic  linear  programs  with  fixed  recourse  (SLPF)  with  discrete  random  elements  have  the  form  in 
(l.l).  They  can  be  formulated  as: 

min  eTx  +Q(x) 

subject  to  Ax  =  b  (1.2) 

z  >  0, 

where 

N 

fi(*)  = 

a*  i 

and  for  each  4*1, ...,  N,  the  recourse  cost  <J(z,  £fc)  is  obtained  by  solving  the  recourse  problem: 

Q(x,  (*)  =  inf  tfy  |  Wy  =  h*  -  I*x,  y  €  }, 

?*  =  (g*,h*,r*), 
pfc  =  probjf  (w)  =  I*). 

Substituting  the  expressions  for  Q  in  (1.2)  ,  we  obtain  a  problem  in  the  form  of  (1.1)  with  W  =  Wk, 
T*  =  «k,  and  pkT*  =  Ck,  for  k  =  1 , ...» AT.  This  problem  has  n  =  no  +  Afni  variables  and  m  =  mo  +  Nm i 
constraints.  The  methods  for  solving  it  include:  the  L-shaped  method,  proposed  by  Van  Slyke  and  Wets  [25); 
the  decomposition  method,  proposed  by  Dantsig  and  Madansky  (11);  and  the  basis  factorisation  method, 
proposed  by  Strasicky  (23),  and  modified  by  Kali  |1?|  and  Wets  (27).  The  first  method  directly  solves  SLPF 
in  form  (1.1),  while  the  other  two  solve  the  dual.  Birge  [3|[4]  discusses  ths  relationship  among  them.  Also 
see  (27|(28|  for  other  references. 


If  Ti  =  =  Tn,  then  (1.1)  has  a  staircase  structure.  This  type  of  problem  is  widely  encountered  in  the 

context  of  dynamical  systems,  i.e.,  discrete  versions  of  continuous  linear  programs  or  linear  control  problems 
[7|[8][13][2l)[25l[28|. 

In  this  paper,  we  study  the  approach  of  Karmarkar’s  projective  algorithm  for  linear  programming  [19] 
to  solve  (l.l).  We  show  that  the  worst-case  performance  bound  can  be  improved  for  this  special  structure 
by  reducing  the  complexity  of  computing  the  projections.  For  large  ~  *io  ~  ni),  this  improvement  is 
an  order  of  magnitude.  We  also  show  that  refinements  of  a  special  case  of  (1.1)  can  be  solved  in  at  most  a 
number  of  iterations  which  is  linear  in  the  additional  data.  These  results  indicate  the  potential  of  interior 
point  methods  in  stochastic  linear  programming,  but  extensive  computational  studies  will  be  required  to 
determine  their  practical  average-case  behavior. 

A  standard  form  version  of  Karmarkar’s  algorithm  is  briefly  described  in  Section  2.  For  L,  the  length  of 
the  input,  this  method  requires  0(nL)  iterations.  The  order  of  each  iteration  is  dominated  by  the  projection 
operation  requiring  generally  0(n3L)  operations.  In  Section  3,  by  means  of  the  Sherman-Morrison- Woodbury 
formula  and  block  matrix  method,  we  explore  the  advantage  of  the  special  structure  of  (1.1)  and  obtain  a 
method  for  the  projection  with  order  <?(AT(nJ  +  njni  +  non3))  =  0(n7L)  if  no  —  nj  ~  N.  A  further 
reduction  is  possible  using  rank-one  updates  of  an  inverse  or  Cholesky  factorisation.  Section  4  presents 
another  variant  of  Karmarkar’s  algorithm  to  solve  the  dual  of  (1.1)  (with  block  angular  structure).  Section 
5  describes  the  implementation  of  the  variants  in  sequential  approximation  methods  and  shows  how  worst 
case  bounds  can  be  improved  in  these  procedures.  Section  6  provides  some  special  cases  where  additional 
computational  savings  are  possible. 

2.  Karmarkar’s  projective  algorithm 

We  briefly  describe  a  standard  form  versionof  Karmarkar’s  projective  algorithm  as  in  Gay(l4],  Ye|30),  or 
Lustig[20j  (with  known  objective  value).  We  choose  the  standard  form  because  we  believe  it  is  more  practical 
for  computation  than  various  canonical  forms.  We  also  assume  an  unknown  objective  value  and  use  Todd 
and  Burrell’s  [24]  method  for  updating  a  lower  bound  on  the  objective  value.  We  use  an  initial  lower  bound 
as  is  often  available  in  practice.  An  alternative  is  Anstreicher’s  [1]  method  to  obtain  an  initial  lower  bound. 


Suppose  ft  linear  program: 


min  cTx 

eubject  to  Ax  =  6,  (2.1) 

x  >  0. 

where  x  €  Rn,  c  €  Zn,  b  6  Zm,  A  €  Zmxn  with  optimal  value  eTx *  =  z*. 

Suppose  that  we  have  a  strictly  interior  feasible  point  a  of  (2.1),  i.e., 

Aa  =  6,  a  >  0,  (2.2) 


a  lower  bound  /°  on  **,  and  that  the  set  of  feasible  solutions  in  (2.1)  is  bounded.  Note  that  if  we  do  not 
have  a  feasible  solution,  then  we  can  solve  a  phase-one  problem  as  in  Karmarkar  [19]. 

Karmarkar’s  algorithm  creates  a  sequence  of  points  z°,  z1, ...,  x*  by  these  steps: 

1.  x°  =  a,k  =  0.  A  lower  bound  1°  <  s*. 

.  * 

£.  If  c7  i*  -  1  is  small  enough,  i.e.,  less  than  a  given  positive  number  c,  then  stop.  Otherwise  go  to  3. 

3.  Let  D  -  diag{xj . *£},  A  :=  [AD,  —6],  and  let  be  the  projection  onto  the  null  space  of  A.  Find 


and 


« -  Pi  (°)  <*  *> 


and  let  /*(/*)  =  min(u,  +  Piti  :»■  1 . »  +  l}. 

If  mV")  <  0,  let  f*+‘  =»  I*. 

Otherwise,  let  fh+1  =  min{u,/tii  :»s>0,i=l,...,»+l}. 

4 ■  Let  cp  =»  u  -  lfc+lu  -  (cTz*  -  l)e/(n  +  1),  where  e  =  (1, . . . ,  1)T  € 

Let 

,  1  a  ep 

»  +  l*  V'n("  +  I)  IM» 

Suppose  /  =  f/jw+i,  where  ?€*",  *,+1  €  R.  Then  **♦*  »  Dj/g**  t.  k  ■  k  +  1,  go  to  1. 
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The  main  computational  effort  at  each  step  k  is  to  compute  the  projections  in  (2.3).  This  effort  is 
generally  0(naL).  Since  this  algorithm  converges  in  0(nL)  steps  (see  Ye[30j),  the  algorithm  has  complexity 
0{n4L3).  Karmarkar  (I9j  uses  a  rank-one  updating  scheme  to  reduce  the  complexity  to  0(n3  hL2).  We 


discuss  this  approach  and  a  related  updating  method  of  Shanno(22)  in  Section  3. 


S.  Computation  Reduction 

We  consider  (1.1)  and  assume  nx  =  =  =  n*  and  mj  =  mj  =  •=  ms  without  loss  of  generality. 

If  an  element  of  a  vector  or  a  matrix  needs  to  be  specified,  we  use  parentheses,  e.g.,  (Ak)iX. 

Let  n  =  no  +  Nnx,m  =  mo  +  NmXt 


-Q'-CMl)' 


and 


A  = 


/  Ao  0 

Ai  Wx 


0  \ 
0 


0  0 
\  As  0  •  WN ) 

Then  (1.1)  can  be  expressed  exactly  the  same  as  (2.1)  for  application  of  Karmarkar’s  algorithm  as  described 
in  Section  2.  Suppose  a  =  (aj , . . . ,  a TN  )r  is  the  current  iteration  point  of  Karmarkar’s  algorithm.  Let 
Di,  =  diag{(ak)i, . . . ,  (afc)„t },  D  =  diag{D0> . . . ,  Ds)-  As  discussed  in  Section  2,  the  main  computational 
work  at  each  step  of  Karmarkar’s  algorithm  is  to  compute  the  projections  in  (2.3).  The  projection  can  be 
written  as 

^»(/-^(A^)-U), 


where  (AAT)  «  ADaj4r  +  bbT  M.  In  this  case,  the  work  is  dominated  by  computing  M~l.  Using 
an  explicit  inverse  would  generally  not  be  the  most  efficient  method  for  computing  the  projection  but  the 
complexity  ia  the  same  in  other  approaches.  In  the  following,  we  use  the  inverse  notation  to  simplify  our 
discussion  of  reduced  complexity  order  for  (1.1).  We  then  describe  how  this  may  be  used  in  practice  to  solve 
systems  with  M. 

Theorem  ».l.  Let  S0  -  -/3  €  5*  -  WhBlWj,  k  *  1 . ff,  5  -  diag{$» . . . ,  S*}.  Then 

S  1  **  diag{5o,51  1,...,5jJ1}.  Let  B%  *  (Asia)  for  fc  «  0,...,N,T5o  *  diag{Do,  1),  h  and  /3  be  unit 


•I  .(  .•  .1  >.|  4»  i.i  .1 


matrices  in  RBo+1  and  Rm°  respectively.  Let 


Gi  =  (Do)~7  +  £  BlS;lBk,G,  =  -BoG^Bl, 


Bo  h 

By  0 

Bn  0 


If  Gx  is  invertible,  then  Gj  and  M  are  invertible  and 


°^)(o  40  U  I)(T 

Proof:  In  Karmarkar’s  algorithm,  a  is  positive,  Le.,  D  is  always  invertible.  Since  Wk  has  full  row  rank,  Sk  is 

invertible  for  jfc  =  1 . N.  Therefore,  S  is  invertible.  Let  D  =  diag{L>0,  fa}-  Let 

,  _  /  BoBo  I2\ 


/  _  _rv  /  BoDo  /a\ 

■-(-i  o°)'D-(  L  ;)-«* 

V  '  \BnDo  Of 


It  is  seen  that  M  —  S  +  UU7.  According  to  the  Sherman- Monrison- Woodbury  formula  [  16],  M  is  invertible 
and 

M~l  =  S~l  -  S~lO(I  +  UT S~ltj)~lUT s~l 

=  s~l  -  S~lU(B~a  +  UTS-lU)~lUTS~l  (3.3) 

=  S'1  -  S~lUG~lUT S'1 
if  and  only  if  (/+  U7 S~lU)  is  invertible,  i.e.,  G  is  invertible. 

Suppose  Gx  is  invertible.  Since  Gx  is  a  symmetric  matrix,  so  is  Gy1,  and  we  can  write  Gj- 1  = 
Gx  1^3GJ’1^3,  where  Gy1*7  is  also  a  symmetric  matrix  (although  possibly  complex).  Since  A0  has  full  row 
rank,  BoGy1^7  also  has  full  row  rank.  Since  no  +  1  >  mo,  Gj  =  — (JJbG71^3)(BoG71^3)r  is  also  invertible. 
It  is  easy  to  see  that 


G(h  or ‘srw/i  o  owor1  o')./"/.  ») 

\°  /,  i,)\o  h)  V°  h) 


■motno-fl' 


Hence,  G  is  invertible.  According  to  (3.5),  M  is  invertible  and  (3.4)  holds.! 

By  using  (3.2),  we  can  now  reduce  the  complexity  of  computing  the  projection  Pjy. 


R 


1 


E 


s 


ft 


« 


Theorem  3.2  (Complexity).  Suppose  G\  is  invertible  then  the  overall  running  time  of  computing  projections 
(2.3)  is  0(JV(nf  +  nfrio  +  n^nji). 

Proof:  At  each  step,  we  have  for  n  =  max{no,  ni), 


work 

calculate  5-1 
calculate  G  i 
calculate  G 2 
calculate  G^1  and  GJ1 
multiply  A  with  a  vector 
multiply  5_1  with  a  vector 
multiply  U  with  a  vector 


number  of  operations 

0(Nn>) 

0(Nnfno  +  ni«§) 

O(ni) 

0(nl) 

0(mn) 

O(mmi) 

0(mr»o) 


Other  efforts  are  small  compared  with  the  operations  above.  According  to  our  assumptions,  the  number 
of  arithmetic  operations  is  dominated  by  0(iV(ni  +  nfno  +  inn^)).  Since  each  arithmetic  operation  needs 
a  precision  of  O(L)  we  obtain  our  conclusion.  ■ 

As  we  stated  above, in  practice,  we  would  not  compute  M~l  explicitly.  The  work  in  (2.3)  is  dominated 
by  the  effort  to  solve  systems  of  the  form 

Mv  =  u,  (3.4) 

using 

v  =  p  -  r,  (3.5) 

where 

Sp  =  u, 

Gt  =  UTp,  (3.6) 

Sr  -  Uq. 

The  systems  in  (3.0)  require  solving  systems  with  5a  .computation  of  Gj  and  Gj  and  solving  systems  with 
G 1  and  Gj.  In  practice,  we  would  find  a  Cholesky  factorisation  of  each  5a,  use  them  to  find  G 1  and  G2  and 
find  Cholesky  factorisations  of  G 1  and  G3. 


Karmarkar  [19]  used  a  rank-one  updating  technique  for  solving  a  canonical  form  linear  program.  In 
this  case,  we  let 

*  =  £i>J+1/*y).  (3-7) 

j=i 

and  update  (Dk) a  only  if 

*a(ft)«/*?+1  *  ( 1/2,2].  (3.8) 

In  our  problem,  for  each  acfc(t)  such  that  (3.8)  holds,  we  would  update  Sk,  S^1,  and  Gv.  In  k  iterations  of  the 
algorithm  (where  lo  =  s*),  we  would  require  0(Jkn°'&)  of  these  updates.  In  this  case,  using  (3.2)  results  in  a 
complexity  of  0((na6n?  +  nn  +  n$)nL3)  for  the  entire  algorithm,  or,  if  N  ~  no  —  Hi,  0(n2  6(?),  compared  to 
Karmarkar’s  general  result  of  0(n3  SLu).  For  unknown  objective  values,  we  could  use  Karmarkar's  canonical 
form  and  (together  with  the  results  in  the  next  section)  obatain  a  similar  reduction. 

Several  practical  problems  occur  with  using  rank-one  updates  of  the  explicit  inverses.  First,  we  cannot 
guarantee  the  lower  bound  estimates  in  our  algorithm  unless  i0  =  **•  We  can,  in  practice,  however,  update 
lk  only  after  several  iterations  or  use  eTxk  (with  different  termination  condition)  in  place  of  lk+l  in  Step  4 
and  still  achieve  low  iteration  counts  (Lustig[20]).  The  second  problem  is  that  we  cannot  guarantee  the  same 
complexity  order  if  we  allow  the  large  values  of  a  used  in  practice  (see,  e.g.,  the  discussion  in  [15]).  Shanno 
[22]  observed  that  the  complexity  can  be  viewed  as  a  function  of  a,  that  no  bound  on  a  appears  necessary, 
and  that  intervals  (1/p,  p]  for  p  >  2  may  also  be  used  in  (3.8)  to  increase  computational  efficiency. 

The  third  problem  in  Karmarkar’s  rank-one  updating  scheme  is  its  use  of  explicit  inverses.  Shanno 
also  Inscribes  how  to  avoid  this  problem  by  using  Fletcher  and  Powell’s  [12]  updating  formula  for  Cholesky 
factorisations  (with  the  same  complexity  order  as  the  inverse  updates).  This  appears  to  achieve  substantial 
increases  in  efficiency.  In  our  case,  we  would  use  this  procedure  to  update  our  Cholesky  factorisations  of  Sk, 
G t,  and  G a,  and,  thereby,  achieve  a  practical  and  efficient  use  of  our  projection  method. 

4.  Dual  Formulation 

The  dual  problem  to  (1.1)  can  also  be  solved  by  a  variant  of  Karmarkar’s  algorithm.  The  program  has 
block-angular  structure  which  leads  to  computational  reductions  similar  to  those  of  Section  3.  An  advantage 


is  also  gained  for  the  use  of  Karmarkar’s  algorithm  in  approximation  schemes. 


The  dual  of  (1.1)  is 


max  +££=!**** 

subject  to  Ajxo  +  I^-i  A***  <  co 

wjxk  <ck,k=l,...,N. 


An  alternative  in  form  (2.1)  can  be  easily  obtained.  Rewrite  (4.1)  as 


min 

-ESLoZ* 

subject  to 

Ylk=o  Ak  yk  =  c0, 

wlyk  =ck,k=  1,. 

•• ,N ; 

i Ik  -  <  0,  k  =  0, . . 

.,W; 

where  aJ  =  (A?, -A?, I), A?  =  (Aj.-Aj, 0),  £  «  (£,-£, 0),  W*  =  (W*  -Wf,/),  k  =  1 . AT,  and 

the  yfc  variables  are  defined  accordingly. 

The  computational  effort  is  again  dominated  by  calculating  M~l  where  M  =  AD2  A?  +  ccT  for  A  the 
constraint  matrix  in  (4.2). 

Note  that  _ 


(ZLoaIdIa*  £d\wx 

ad2a=  "ta*  o 

-r  =  -  0 

V  W^A*  0 


atnd%wn 


WTNDj,WN 


Now  observe  that  M  =  P  +  C/VT  where 

/  n2T 


/Ei»o AkD2hAk  0 
p\  WlDlM  WTlD2Wl  0 


Kd2nan 


WT„D2„Wn 


/  c0 

0  ci 
0  c/f 


0  _  co 

WiD\M  a 

FT?!"  nQ  "V 


\WlOjfAff  cs 

The  Sherman- Morrison- Woodbury  formula  yields 


A/~‘  =  P~l  -  P~lU{I  +  VT P~lU)-lVT P~l . 
9 


I 


1 

K 

I 

V 

?' 

y 


i 


I 

H 


H  =  ^ATkD\Ak, 


Hk  =  W*D\Ak,k  = 


Jk  =  wZDlWk,k=l,...,N, 


P-1  = 


H~l 

-jzihxh~i 


(I+VTP~1U)  =  (  I  2*=i  £*Jk  1  Yhc=i(HkJk  lHkH  1co  —  HkJklc0)  \ 

1  1  \cT0H-'-Y.LlcTkJ;lBkH-'  I  +  cTH-'c0  +  Z,?=1(clJ;1HkH-'co  +  cTj^ck)) 


Equations  (4.3),  (4.4),  and  (4.5)  can  then  be  used  to  develop  another  projection  method  for  (4.2). 
Theorem  4.2.  The  overall  running  time  of  the  projection  using  (4.3-5)  is  0(JV(nf  +  nfn<,  +  nir»§)). 


Proof:  The  computational  effort  in  each  step  is: 
work 


Calculate  H~x 


Calculate  Jkl,k  =  1 . . . , N 
Calculate  HkH~x,  k  «  1, . . . ,  N 
Calculate  H%Jkl,k  =  1 ,..., JV 
Calculate  (/  +  VTP~lU)~l 
Multiply  A  with  a  vector 
Multiply  P~ 1  with  a  vector 
Multiply  VT  with  a  vector 


number  of  operations 
O(JVmng) 

0[Nn\) 

O(Nninl) 

O(Nnlno) 

0(n§) 

0(mno  +  TVn^mi  +  nj)) 
O(nn) 

O(non) 


AU  other  effort  is  dominated  by  these  operations.  The  number  of  operations  is  0(N(n\  +  njno  +  rijn*)), 


and  the  result  is  obtained  as  in  Theorem  3.2.1 


Theorem  4.2  gives  the  same  complexity  as  Theorem  3.2,  so,  again  for  N  ~  no  ~  ni,  we  have  0(n3L2) 
complexity  for  the  full  algorithm.  We  can  again  use  Cholesky  factorizations  to  make  the  projection  with 
(4.3-5)  practical. 

5.  Dual  Use  in  Approximations 

The  dual  provides  an  additional  benefit.  In  using  (3.1)  to  approximate  a  stochastic  linear  program  with 
more  than  N  realizations  of  the  random  vector  bk,  additional  realizations  fejv-t- i ,  •  •  • ,  6jv+i  are  added  to  (3.1) 
and  the  probabilities  p*,  k  —  1, . . . ,  N  + 1  are  updated.  (See  Birge  and  Wets  [6]  for  procedures  to  do  this.)  A 
feasible  solution  to  the  new  problem  (3.1)  with  N'  =*  N  + 1  is  not  readily  found  but  given  a  feasible  solution 
(4.2)  another  can  be  found  easily  if  T1  =  T2  =  •  •  •  =  TN  =  T  and  q1  =  q2  =  •  •  •  =  qN  =  q.  Under  these 
assumptions  let  the  old  solution  be  (yo>d,  ■  ■  • ,  y$d)  where 

WTylld  =pkq,k=l,...,N, 

and 

N 

A? ..old  .  \  <7 5^"  old  _ 

Ao  Vo  +  2^>  1  y*  ~  C°> 

*=1 

where  T  =  (T,  —  T,0). 

Define  a  new  solution  by 

,,uvm  _  old 

Vo  —  Vo  i 

N  J5 

vnew  =  (Etfld)prr.  *  =  !,•••.*+<• 

»=1 

The  solution  in  (5.1)  is  then  feasible  in  (4.2)  with  objective  value  6g  ygld  +  Y!>k=i  ^k  (Sili  y»?ld)P*'w- 
The  total  effort  for  solving  the  new  (4.2)  should  be  less  than  resolving  (3.1)  because  no  procedures  to  find 
a  feasible  interior  point  an  necessary  (assuming,  of  course,  that  yold  is  interior).  The  objective  is  generally 
close  to  the  old  optimal  value  of  (4.2)  and  should  be  close  to  the  new  optimal  value.  The  effort  factor  for 
the  number  of  iterations  to  solve  the  new  (4.2)  should  therefore  be  much  smaller  than  O(nL). 

A  similar  procedure  can  be  used  in  the  primal  problem  (3.1)  if  the  randomness  is  limited  to  the  objective 
function,  i.e.,  T1  =  T3  =  ■  ■  ■  =  TN  =  T  and  h1  =  ...  =  hN  =  h  but  q *  /  q>  for  all »'  j.  Again, a  starting 


feasible  solution  can  be  easily  found  by  letting 


_new 

xk 


x%u,k  =  0,...,N, 


x“ew  =  argminIoid{gJzfld)t  =  =  N  +  1, . . . ,  N  +  l. 

Now  let  p  =  p“ew,  2oid  be  the  objective  of  the  old  (3.1)  with  z£ld  and  let 


(5.2) 


ff+t 


»0  _  jr.T- new  .  \  '  _ne»J_ new 

*new  “PC  *0  +  2-e  P1  ^ Xi  ’ 

j=N+ 1 


N 


The  new  initial  objective  value  in  (3.1)  is  then 

Znc*  =  (1  -  p)cTx Sld  +  Y^P*tmqixfd  +  z°tw, 

and  we  wish  to  find 

=  min{(l  -  p)crzgld  +  ]Tp“eV*y  +  z'N), 


)=i 


N 


(5.3) 


j=i 


an 


where  z'N  =  pcTx0  +  p"*w^zy.  Let  2*  be  the  optimal  value  of  the  old  (3.1)  and  let  xold  be 

interior  solution  of  the  old  (3.1)  such  that  20ld  <**+€.  Now  assume  that  p“*w/(  1  -  p),  i.e.  that  the 
old  probabilities  have  the  same  relative  values  in  the  new  problem,  so  that  (1  -  p)z*  is  the  minimum  of 
(1  -  p)cT*0  +  p“*w?JXj.  We  then  have 

*new  £  (1  -  P)*o!d  +  max  g'N 

(5.4) 


*new  t  (i  -  p)(*oid  -  c)  +  min  z'N. 

Let  L'  be  the  si se  of  the  data  for  the  objective  and  constraints  relating  only  to  z'N.  This  yields  the 


following  result. 


Theorem  5.1  The  overall  running  time  of  Algorithm  3.1  starting  from  a  solution  defined  by  (5.2)  for  N  + 1 
realisations  of  the  objective  g*,  assuming  T1  =  =  TN+l  =  T,  h1  =  -  nN+l  =  h,  p"*w/(l  -  p)  =  p?ld, 

and  «  <  2^,  is  0[(N  +  l)(nf  +  nfno  +  ning)nLL'). 

Proof:  The  number  of  operations  per  iteration  follows  from  Theorem  3.2.  The  number  of  iterations  0(nL') 
follows  from  (5.4)  for  2°^’*  >  z'N  >  —2°^). 


Theorem  5.1  provides  a  worst-case  bound  that  can  be  used  in  solving  successively  larger  stochastic 
linear  programs.  The  entire  approximating  process  can  then  be  seen  as  a  variable  dimension  variant  of 


Karmarkar’s  algorithm  in  which  the  number  of  iterations  is  0(lnL')  where  l  additional  problems  of  size  L' 
are  added  after  the  initial  solution  of  (3.1). 

This  result  contrasts  with  the  alternative  (simplex-based)  approaches  to  stochastic  linear  programs,  for 
which,  the  complexity  of  solving  each  new  problem  depends  on  the  size  of  the  entire  problem.  Good  bounds 
on  t'n  that  would  often  be  available  in  practice  could  further  reduce  the  computational  effort. 

For  a  single  SLPF,  we  can  also  compare  the  projective  algorithm  with  simplex-  based  methods  if  we 
assume  average  cases  for  iteration  counts.  The  best  available  results  for  stochastic  linear  programming 
appear  to  be  in  uses  of  the  L-shaped  method  of  Van  Slyke  and  Wets  (Birge  [2]).  It  appears  that  the 
running  time  is  typically  then  0(ATn3)  =  0(n2)  (Birge  and  Louveaux  [5]).  If  we  use  rank-one  updates  in  the 
projective  method  and  assume  an  average  O(logn)  iterations  as  has  been  conjectured  (Karmarkar  [18]),  then 
the  projective  algorithm  has  running  time  O(ni  slogn)  (ignoring  L  which  is  assumed  in  both  estimates). 
This  observation  combined  with  the  potential  for  sequential  approximation  indicates  the  promise  for  the 
projective  algorithm  in  stochastic  linear  programming.  Further  computational  experience  will  be  required 
to  bear  this  out. 


6.  Special  Cases 

The  work  of  the  variants  of  Karmarkar’s  algorithm  is  reduced  for  the  following  special  cases. 

(1)  Simple  recourse  [28) 

In  this  case,  W  =  (1,-1),  where  I  is  the  unit  matrix  in  and  nj  =  2 Write  **  = 

(xk  >xk)>xk  >xh  e  Rm,i  for  *  -  !»•••.  N  in  (3.1).  Let  Dk  =  diag{I)£,i?fc}  in  (3.2),  where  £>£,£>*  € 
Smi  Xm* ,  for  k  =  1, . . . ,  N.  Then  in  Theorem  3.1,  we  have 


5*  =  WDlWT  =  (Z>+)3  +  (DZ)\ 

Sfc-1  =  diag{((«+)»  +  (a,-)?)"1 . ((°k  )m,  +  (®a  )mi)-1}> 


(6.1) 


where  Df  =  diag^a^, 


•>(°fc)m,}  nnd  Dl  =  diag{(o*  )l( . . . ,  (oj  )„,},  for  k  =  1  Therefore, 


S'1  can  be  easily  calculated  in  (3.2). 


The  work  is  reduced  but  the  order  of  the  algorithm  is  the  tame.  Similar  modifications  can  be  made  in 


the  dual  approach. 


(S)  Network  recourse  (26) 
In  this  case, 


where  E  is  the  network  node-arc  incidence  matrix.  The  calculation  work  of  5s  and  1  can  still  be  reduced 
but  again  the  order  is  the  same. 
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